## subject & microphone list ----
## fomri arrive (range):
fomri.arrive.lb <- "2018-05-01"
fomri.arrive.ub <- "2018-05-31"
## definite fomri
fomri.switch <- "2018-06-28"
## exceptions:
## DMCC515268_baseline DMCC3 (June 30th, 2018) and DMCC6904377_reactive DMCC2 (July 1st, 2018)
dates1.is.microoptics <- dates1[c("BAS.day", "PRO.day", "REA.day")] < fomri.arrive.lb
dates1.is.fomri <- dates1[c("BAS.day", "PRO.day", "REA.day")] > fomri.switch
all(dates1.is.microoptics, na.rm = TRUE) ## all microoptics## [1] TRUEdates2.is.microoptics <- dates2[c("session1", "session2", "session3")] < fomri.arrive.lb
dates2.is.fomri <- dates2[c("session1", "session2", "session3")] > fomri.switch
dates2$all.microoptics <- apply(dates2.is.microoptics, 1, function(x) all(x))
dates2$all.fomri <- apply(dates2.is.fomri, 1, function(x) all(x))
subj.microoptics <- dates2 %>% filter(all.microoptics & subj != "DMCC6904377") %>% .$subj
subj.microoptics <- union(subj.microoptics, dates1$subj)
subj.fomri <- dates2 %>% filter(all.fomri & subj != "DMCC6904377") %>% .$subj
stroop$mic <- ifelse(
stroop$subj %in% subj.fomri,
"fomri",
ifelse(
stroop$subj %in% subj.microoptics,
"micro",
"unknown"
)
)
table(stroop$subj, stroop$mic)##
## fomri micro unknown
## 107321 0 0 672
## 115825 0 672 0
## 123117 0 672 0
## 130114 0 672 0
## 130518 0 672 0
## 132017 0 672 0
## 138837 0 672 0
## 141422 0 672 0
## 150423 0 672 0
## 158136 0 672 0
## 160830 0 0 672
## 161832 0 0 672
## 165032 0 672 0
## 171330 0 456 0
## 178243 0 645 0
## 178647 0 672 0
## 178950 0 672 0
## 182436 0 0 672
## 197449 0 672 0
## 203418 0 672 0
## 204319 0 672 0
## 233326 0 672 0
## 300618 0 672 0
## 317332 0 672 0
## 346945 0 672 0
## 352738 0 432 0
## 393550 0 672 0
## 448347 672 0 0
## 580650 0 672 0
## 594156 0 672 0
## 601127 0 672 0
## 765864 0 672 0
## 814649 0 672 0
## 843151 0 456 0
## 849971 0 672 0
## 873968 672 0 0
## 877168 0 672 0
## DMCC1328342 672 0 0
## DMCC1596165 672 0 0
## DMCC1971064 0 672 0
## DMCC2442951 672 0 0
## DMCC2803654 0 672 0
## DMCC2834766 0 672 0
## DMCC3062542 672 0 0
## DMCC4191255 0 672 0
## DMCC4260551 0 0 672
## DMCC4854075 0 672 0
## DMCC5009144 456 0 0
## DMCC5195268 0 672 0
## DMCC5775387 0 672 0
## DMCC6371570 672 0 0
## DMCC6418065 0 672 0
## DMCC6484785 0 0 672
## DMCC6661074 0 432 0
## DMCC6671683 0 672 0
## DMCC6705371 0 672 0
## DMCC6721369 0 672 0
## DMCC6904377 0 0 456
## DMCC7297690 432 0 0
## DMCC7921988 0 0 672
## DMCC8033964 0 456 0
## DMCC8050964 672 0 0
## DMCC8078683 0 0 672
## DMCC9441378 672 0 0
## DMCC9478705 0 672 0
## DMCC9953810 672 0 0stroop %>%
ggplot(aes(subj, rt, color = mic)) +
facet_grid(vars(session)) +
geom_point(position = position_jitter(width = 0.1), size = 0.5, alpha = 0.3) +
geom_boxplot(position = position_nudge(1/3), notch = TRUE, width = 0.2) +
scale_color_brewer(type = "qual", palette = 2) +
theme(axis.text.x = element_text(angle = 90, hjust = 0, vjust = 0))## Warning: Removed 422 rows containing non-finite values (stat_boxplot).## Warning: Removed 422 rows containing missing values (geom_point).stroop %>%
group_by(mic, subj) %>%
summarize(freq = sum(rt == 0, na.rm = TRUE)) %>%
ggplot(aes(mic, freq)) +
geom_point(size = 2) +
labs(y = "frequency of rt == 0 per subj")stroop %>%
group_by(mic, subj, error) %>%
summarize(freq = sum(error, na.rm = TRUE)) %>%
ggplot(aes(mic, freq)) +
geom_point(size = 2, position = position_jitter(width = 0.1), alpha = 0.4) +
labs(y = "frequency of errors per subj")stroop %>%
group_by(mic, subj, acc.final) %>%
summarize(freq = n()) %>%
ggplot(aes(mic, freq)) +
facet_grid(cols = vars(acc.final)) +
geom_point(size = 2, position = position_jitter(width = 0.1), alpha = 0.4) +
labs(y = "frequency of acc.final codings per subj")stroop %>%
ggplot(aes(sample = rt, color = mic)) +
stat_qq(alpha = 0.3, size = 0.15) +
stat_qq_line(size = 0.25) +
facet_wrap(vars(subj)) +
scale_color_brewer(type = "qual", palette = 2)## Warning: Removed 422 rows containing non-finite values (stat_qq).## Warning: Removed 422 rows containing non-finite values (stat_qq_line).stroop.mean.sd <- stroop %>%
group_by(mic, subj, session) %>%
summarize(rt.mean = mean(rt, na.rm = TRUE), rt.sd = sd(rt, na.rm = TRUE))
grid.arrange(
stroop.mean.sd %>%
ggplot(aes(mic, rt.mean)) +
facet_grid(vars(session)) +
geom_point(position = position_jitter(width = 0.1)) +
geom_boxplot(position = position_nudge(1/3), notch = FALSE, width = 0.2),
stroop.mean.sd %>%
ggplot(aes(mic, rt.sd)) +
facet_grid(vars(session)) +
geom_point(position = position_jitter(width = 0.1)) +
geom_boxplot(position = position_nudge(1/3), notch = FALSE, width = 0.2)
)stroop.rt <- stroop %>% filter(!is.na(rt), rt > 200, acc == 1)
# g <- paste0(stroop.rt$subj, "_", stroop.rt$session)
# l <- split(stroop.rt, g)
#
# fits <- lapply(l, function(d) lm(rt ~ trial.type, d))
# stroop.rt$resids <- unlist(lapply(fits, resid), use.names = FALSE)
# bs <- do.call(rbind, lapply(fits, coef))
# bs <- cbind(bs, reshape2::colsplit(rownames(bs), "_", c("subj", "session"))
stroop.rt %<>%
mutate(
bas = as.numeric(session == "bas"),
pro = as.numeric(session == "pro"),
rea = as.numeric(session == "rea"),
congr = as.numeric(trial.type == "c"),
incon = as.numeric(trial.type == "i"),
micro = as.numeric(mic == "micro"),
fomri = as.numeric(mic == "fomri")
)
fit.noint <- lmer(
rt ~
0 +
congr:bas:micro + congr:pro:micro + congr:rea:micro +
incon:bas:micro + incon:pro:micro + incon:rea:micro +
congr:bas:fomri + congr:pro:fomri + congr:rea:fomri +
incon:bas:fomri + incon:pro:fomri + incon:rea:fomri +
(
0 +
congr:bas + congr:pro + congr:rea +
incon:bas + incon:pro + incon:rea
| subj
),
stroop.rt,
subset = mic != "unknown",
control = lmerControl(optimizer = "bobyqa")
)
# fits.noint <- allFit(fit.noint)
# summary(fits.noint)
# fit.noint <- fits.noint$bobyqa
summary(fit.noint)## Linear mixed model fit by REML ['lmerMod']
## Formula: rt ~ 0 + congr:bas:micro + congr:pro:micro + congr:rea:micro +
## incon:bas:micro + incon:pro:micro + incon:rea:micro + congr:bas:fomri +
## congr:pro:fomri + congr:rea:fomri + incon:bas:fomri + incon:pro:fomri +
## incon:rea:fomri + (0 + congr:bas + congr:pro + congr:rea +
## incon:bas + incon:pro + incon:rea | subj)
## Data: stroop.rt
## Control: lmerControl(optimizer = "bobyqa")
## Subset: mic != "unknown"
##
## REML criterion at convergence: 467003
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -5.019 -0.532 -0.123 0.347 12.417
##
## Random effects:
## Groups Name Variance Std.Dev. Corr
## subj congr:bas 14727 121
## congr:pro 12694 113 0.71
## congr:rea 18702 137 0.71 0.84
## bas:incon 19455 139 0.89 0.66 0.60
## pro:incon 16440 128 0.72 0.97 0.79 0.76
## rea:incon 20603 144 0.75 0.87 0.94 0.75 0.90
## Residual 27957 167
## Number of obs: 35655, groups: subj, 57
##
## Fixed effects:
## Estimate Std. Error t value
## congr:bas:micro 830.0 18.5 44.9
## congr:micro:pro 821.6 17.1 48.2
## congr:micro:rea 810.7 20.6 39.3
## bas:micro:incon 962.4 21.3 45.2
## micro:pro:incon 900.2 19.2 46.8
## micro:rea:incon 909.4 21.6 42.0
## congr:bas:fomri 734.6 36.0 20.4
## congr:pro:fomri 718.8 33.0 21.8
## congr:rea:fomri 767.1 40.2 19.1
## bas:incon:fomri 860.2 41.3 20.8
## pro:incon:fomri 810.0 37.2 21.8
## rea:incon:fomri 847.6 42.1 20.1
##
## Correlation of Fixed Effects:
## cngr:bs:m cngr:mcr:p cngr:mcr:r bs:mc: mcr:p: mcr:r: cngr:bs:f
## cngr:mcr:pr 0.684
## congr:mcr:r 0.689 0.813
## bas:mcr:ncn 0.878 0.635 0.579
## micr:pr:ncn 0.698 0.945 0.774 0.737
## micro:r:ncn 0.724 0.851 0.923 0.728 0.883
## cngr:bs:fmr 0.000 0.000 0.000 0.000 0.000 0.000
## cngr:pr:fmr 0.000 0.000 0.000 0.000 0.000 0.000 0.682
## congr:r:fmr 0.000 0.000 0.000 0.000 0.000 0.000 0.681
## bas:ncn:fmr 0.000 0.000 0.000 0.000 0.000 0.000 0.878
## pro:ncn:fmr 0.000 0.000 0.000 0.000 0.000 0.000 0.696
## rea:ncn:fmr 0.000 0.000 0.000 0.000 0.000 0.000 0.718
## cngr:p: cngr:r: bs:nc: pr:nc:
## cngr:mcr:pr
## congr:mcr:r
## bas:mcr:ncn
## micr:pr:ncn
## micro:r:ncn
## cngr:bs:fmr
## cngr:pr:fmr
## congr:r:fmr 0.808
## bas:ncn:fmr 0.634 0.574
## pro:ncn:fmr 0.946 0.769 0.734
## rea:ncn:fmr 0.848 0.922 0.723 0.879bnames <- rownames(coef(summary(fit.noint)))
m <- data.frame(
bas = grepl("bas", bnames),
pro = grepl("pro", bnames),
rea = grepl("rea", bnames),
congr = grepl("congr", bnames),
incon = grepl("incon", bnames),
fomri = grepl("fomri", bnames),
micro = grepl("micro", bnames)
)
rownames(m) <- bnames
contrast <- rbind(
"fomri-micro|bas" = (m$fomri - m$micro) * m$bas / 2,
"fomri-micro|pro" = (m$fomri - m$micro) * m$pro / 2,
"fomri-micro|rea" = (m$fomri - m$micro) * m$rea / 2,
"stroop:(fomri-micro)|bas" = (m$incon - m$congr) * (m$fomri - m$micro) * m$bas / 2,
"stroop:(fomri-micro)|pro" = (m$incon - m$congr) * (m$fomri - m$micro) * m$pro / 2,
"stroop:(fomri-micro)|rea" = (m$incon - m$congr) * (m$fomri - m$micro) * m$rea / 2
)
colnames(contrast) <- bnames
results <- glht(fit.noint, contrast, alternative = "two.sided")
summary(results, adjusted(type = "none"))##
## Simultaneous Tests for General Linear Hypotheses
##
## Fit: lmer(formula = rt ~ 0 + congr:bas:micro + congr:pro:micro + congr:rea:micro +
## incon:bas:micro + incon:pro:micro + incon:rea:micro + congr:bas:fomri +
## congr:pro:fomri + congr:rea:fomri + incon:bas:fomri + incon:pro:fomri +
## incon:rea:fomri + (0 + congr:bas + congr:pro + congr:rea +
## incon:bas + incon:pro + incon:rea | subj), data = stroop.rt,
## control = lmerControl(optimizer = "bobyqa"), subset = mic !=
## "unknown")
##
## Linear Hypotheses:
## Estimate Std. Error z value Pr(>|z|)
## fomri-micro|bas == 0 -98.83 42.12 -2.35 0.019 *
## fomri-micro|pro == 0 -96.56 39.00 -2.48 0.013 *
## fomri-micro|rea == 0 -52.71 45.35 -1.16 0.245
## stroop:(fomri-micro)|bas == 0 -3.36 11.12 -0.30 0.762
## stroop:(fomri-micro)|pro == 0 6.27 6.93 0.91 0.365
## stroop:(fomri-micro)|rea == 0 -9.08 9.18 -0.99 0.322
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## (Adjusted p values reported -- none method)Summary
cl1 <- lmeControl(
maxIter = 100000, msMaxIter = 100000, niterEM = 100000,
msMaxEval = 100000, tolerance = 0.000001, msTol = 0.0000001, returnObject = TRUE,
minAbsParApVar = 0.05, opt = c("nlminb"), optimMethod = "BFGS"
)
## proactive
m.hom.pro <- lme(
rt ~ trial.type * mic,
random = ~ 1 | subj,
data = stroop.rt %>% filter(mic != "unknown"),
subset = session == "pro",
control = cl1
)
m.het.mic.pro <- update(m.hom.pro, . ~ ., weights = varIdent(form = ∼ 1 | mic))
m.het.subj.pro <- update(m.hom.pro, . ~ ., weights = varIdent(form = ∼ 1 | subj))
summary(m.hom.pro)## Linear mixed-effects model fit by REML
## Data: stroop.rt %>% filter(mic != "unknown")
## Subset: session == "pro"
## AIC BIC logLik
## 157167 157211 -78578
##
## Random effects:
## Formula: ~1 | subj
## (Intercept) Residual
## StdDev: 122 169.6
##
## Fixed effects: rt ~ trial.type * mic
## Value Std.Error DF t-value p-value
## (Intercept) 718.8 35.69 11915 20.140 0.0000
## trial.typei 91.1 7.09 11915 12.839 0.0000
## micmicro 102.7 40.17 55 2.557 0.0134
## trial.typei:micmicro -12.4 8.00 11915 -1.546 0.1220
## Correlation:
## (Intr) trl.ty micmcr
## trial.typei -0.132
## micmicro -0.888 0.117
## trial.typei:micmicro 0.117 -0.887 -0.132
##
## Standardized Within-Group Residuals:
## Min Q1 Med Q3 Max
## -4.1408 -0.5416 -0.1393 0.3415 12.1557
##
## Number of Observations: 11974
## Number of Groups: 57summary(m.het.mic.pro)## Linear mixed-effects model fit by REML
## Data: stroop.rt %>% filter(mic != "unknown")
## Subset: session == "pro"
## AIC BIC logLik
## 157121 157173 -78553
##
## Random effects:
## Formula: ~1 | subj
## (Intercept) Residual
## StdDev: 122 173.4
##
## Variance function:
## Structure: Different standard deviations per stratum
## Formula: ~1 | mic
## Parameter estimates:
## micro fomri
## 1.0000 0.8939
## Fixed effects: rt ~ trial.type + mic + trial.type:mic
## Value Std.Error DF t-value p-value
## (Intercept) 718.8 35.61 11915 20.183 0.0000
## trial.typei 91.1 6.48 11915 14.052 0.0000
## micmicro 102.7 40.11 55 2.561 0.0132
## trial.typei:micmicro -12.4 7.50 11915 -1.649 0.0992
## Correlation:
## (Intr) trl.ty micmcr
## trial.typei -0.121
## micmicro -0.888 0.107
## trial.typei:micmicro 0.104 -0.864 -0.124
##
## Standardized Within-Group Residuals:
## Min Q1 Med Q3 Max
## -4.0508 -0.5440 -0.1399 0.3416 11.8919
##
## Number of Observations: 11974
## Number of Groups: 57summary(m.het.subj.pro)## Linear mixed-effects model fit by REML
## Data: stroop.rt %>% filter(mic != "unknown")
## Subset: session == "pro"
## AIC BIC logLik
## 153344 153802 -76610
##
## Random effects:
## Formula: ~1 | subj
## (Intercept) Residual
## StdDev: 121.8 159.9
##
## Variance function:
## Structure: Different standard deviations per stratum
## Formula: ~1 | subj
## Parameter estimates:
## 115825 123117 130114 130518 132017 138837
## 1.0000 1.3684 0.9483 0.7232 1.1539 1.0273
## 141422 150423 158136 165032 171330 178243
## 1.0303 1.3350 0.6193 1.6752 0.9424 1.5864
## 178647 178950 197449 203418 204319 233326
## 1.1050 0.7494 1.1987 0.7306 1.8857 2.0556
## 300618 317332 346945 352738 393550 448347
## 0.5359 0.6395 0.8655 1.6850 0.7852 1.1969
## 580650 594156 601127 765864 814649 843151
## 2.0928 0.6314 1.1921 0.5608 0.9475 1.1984
## 849971 873968 877168 DMCC1328342 DMCC1596165 DMCC1971064
## 1.8200 0.7821 1.3117 0.4071 1.6184 0.8052
## DMCC2442951 DMCC2803654 DMCC2834766 DMCC3062542 DMCC4191255 DMCC4854075
## 0.4382 0.6706 0.7042 0.9748 0.6480 0.6854
## DMCC5009144 DMCC5195268 DMCC5775387 DMCC6371570 DMCC6418065 DMCC6661074
## 0.5089 0.5855 1.3942 0.9953 0.5380 0.7065
## DMCC6671683 DMCC6705371 DMCC6721369 DMCC7297690 DMCC8033964 DMCC8050964
## 0.5219 0.7398 0.5582 0.7739 0.5176 0.8544
## DMCC9441378 DMCC9478705 DMCC9953810
## 1.4022 0.8540 0.9206
## Fixed effects: rt ~ trial.type + mic + trial.type:mic
## Value Std.Error DF t-value p-value
## (Intercept) 735.5 35.44 11915 20.756 0.0000
## trial.typei 65.2 4.67 11915 13.940 0.0000
## micmicro 88.6 39.90 55 2.221 0.0305
## trial.typei:micmicro 8.7 5.45 11915 1.589 0.1121
## Correlation:
## (Intr) trl.ty micmcr
## trial.typei -0.087
## micmicro -0.888 0.078
## trial.typei:micmicro 0.075 -0.858 -0.091
##
## Standardized Within-Group Residuals:
## Min Q1 Med Q3 Max
## -4.9122 -0.6275 -0.1734 0.4266 10.1056
##
## Number of Observations: 11974
## Number of Groups: 57anova(m.hom.pro, m.het.mic.pro, m.het.subj.pro)## Model df AIC BIC logLik Test L.Ratio p-value
## m.hom.pro 1 6 157167 157211 -78578
## m.het.mic.pro 2 7 157121 157173 -78553 1 vs 2 48 <.0001
## m.het.subj.pro 3 62 153344 153802 -76610 2 vs 3 3887 <.0001## baseline
m.hom.bas <- lme(
rt ~ trial.type * mic,
random = ~ 1 | subj,
data = stroop.rt %>% filter(mic != "unknown"),
subset = session == "bas",
control = cl1
)
m.het.mic.bas <- update(m.hom.bas, . ~ ., weights = varIdent(form = ∼ 1 | mic))
m.het.subj.bas <- update(m.hom.bas, . ~ ., weights = varIdent(form = ∼ 1 | subj))
summary(m.hom.bas)## Linear mixed-effects model fit by REML
## Data: stroop.rt %>% filter(mic != "unknown")
## Subset: session == "bas"
## AIC BIC logLik
## 144451 144495 -72220
##
## Random effects:
## Formula: ~1 | subj
## (Intercept) Residual
## StdDev: 123.1 155.2
##
## Fixed effects: rt ~ trial.type * mic
## Value Std.Error DF t-value p-value
## (Intercept) 742.0 37.31 11101 19.889 0.0000
## trial.typei 130.4 6.80 11101 19.186 0.0000
## micmicro 86.1 41.91 51 2.054 0.0451
## trial.typei:micmicro 3.1 7.66 11101 0.409 0.6829
## Correlation:
## (Intr) trl.ty micmcr
## trial.typei -0.060
## micmicro -0.890 0.053
## trial.typei:micmicro 0.053 -0.888 -0.060
##
## Standardized Within-Group Residuals:
## Min Q1 Med Q3 Max
## -4.9472 -0.5558 -0.1074 0.4048 11.8554
##
## Number of Observations: 11156
## Number of Groups: 53summary(m.het.mic.bas)## Linear mixed-effects model fit by REML
## Data: stroop.rt %>% filter(mic != "unknown")
## Subset: session == "bas"
## AIC BIC logLik
## 144351 144402 -72168
##
## Random effects:
## Formula: ~1 | subj
## (Intercept) Residual
## StdDev: 123.1 160.2
##
## Variance function:
## Structure: Different standard deviations per stratum
## Formula: ~1 | mic
## Parameter estimates:
## micro fomri
## 1.0000 0.8417
## Fixed effects: rt ~ trial.type + mic + trial.type:mic
## Value Std.Error DF t-value p-value
## (Intercept) 742.0 37.26 11101 19.917 0.0000
## trial.typei 130.4 5.91 11101 22.083 0.0000
## micmicro 86.1 41.87 51 2.056 0.0449
## trial.typei:micmicro 3.1 6.94 11101 0.451 0.6521
## Correlation:
## (Intr) trl.ty micmcr
## trial.typei -0.052
## micmicro -0.890 0.047
## trial.typei:micmicro 0.044 -0.851 -0.054
##
## Standardized Within-Group Residuals:
## Min Q1 Med Q3 Max
## -4.7921 -0.5606 -0.1103 0.4044 11.4845
##
## Number of Observations: 11156
## Number of Groups: 53summary(m.het.subj.bas)## Linear mixed-effects model fit by REML
## Data: stroop.rt %>% filter(mic != "unknown")
## Subset: session == "bas"
## AIC BIC logLik
## 142054 142478 -70969
##
## Random effects:
## Formula: ~1 | subj
## (Intercept) Residual
## StdDev: 122.4 180.7
##
## Variance function:
## Structure: Different standard deviations per stratum
## Formula: ~1 | subj
## Parameter estimates:
## 115825 123117 130114 130518 132017 138837
## 1.0000 1.1149 0.9185 0.5873 0.8056 1.1128
## 141422 150423 158136 165032 178243 178647
## 0.7795 0.8202 0.6938 0.6761 1.5863 1.1238
## 178950 197449 203418 204319 233326 300618
## 0.7375 1.7405 0.8050 0.8296 1.5524 0.6524
## 317332 346945 352738 393550 448347 580650
## 0.5346 1.0779 0.6380 0.7650 1.1714 0.6154
## 594156 601127 765864 814649 849971 873968
## 0.5954 1.1536 0.5274 1.0802 1.3437 0.7165
## 877168 DMCC1328342 DMCC1596165 DMCC1971064 DMCC2442951 DMCC2803654
## 1.1406 0.4654 0.6111 0.8529 0.3401 0.7688
## DMCC2834766 DMCC3062542 DMCC4191255 DMCC4854075 DMCC5195268 DMCC5775387
## 0.6989 0.6309 0.5408 0.3850 0.7070 0.7052
## DMCC6371570 DMCC6418065 DMCC6661074 DMCC6671683 DMCC6705371 DMCC6721369
## 0.7652 0.7075 0.7614 0.5515 0.5280 0.6040
## DMCC7297690 DMCC8050964 DMCC9441378 DMCC9478705 DMCC9953810
## 0.6113 0.9713 0.9493 0.8234 0.6086
## Fixed effects: rt ~ trial.type + mic + trial.type:mic
## Value Std.Error DF t-value p-value
## (Intercept) 747.9 37.03 11101 20.196 0.0000
## trial.typei 112.3 4.76 11101 23.604 0.0000
## micmicro 82.6 41.61 51 1.985 0.0526
## trial.typei:micmicro 12.2 5.60 11101 2.175 0.0297
## Correlation:
## (Intr) trl.ty micmcr
## trial.typei -0.042
## micmicro -0.890 0.038
## trial.typei:micmicro 0.036 -0.849 -0.044
##
## Standardized Within-Group Residuals:
## Min Q1 Med Q3 Max
## -4.4068 -0.6214 -0.1296 0.4707 10.5597
##
## Number of Observations: 11156
## Number of Groups: 53anova(m.hom.bas, m.het.mic.bas, m.het.subj.bas)## Model df AIC BIC logLik Test L.Ratio p-value
## m.hom.bas 1 6 144451 144495 -72220
## m.het.mic.bas 2 7 144351 144402 -72168 1 vs 2 102.7 <.0001
## m.het.subj.bas 3 58 142054 142478 -70969 2 vs 3 2398.8 <.0001## reactive
m.hom.rea <- lme(
rt ~ trial.type * mic,
random = ~ 1 | subj,
data = stroop.rt %>% filter(mic != "unknown"),
subset = session == "rea",
control = cl1
)
m.het.mic.rea <- update(m.hom.rea, . ~ ., weights = varIdent(form = ∼ 1 | mic))
m.het.subj.rea <- update(m.hom.rea, . ~ ., weights = varIdent(form = ∼ 1 | subj))
summary(m.hom.rea)## Linear mixed-effects model fit by REML
## Data: stroop.rt %>% filter(mic != "unknown")
## Subset: session == "rea"
## AIC BIC logLik
## 165795 165839 -82891
##
## Random effects:
## Formula: ~1 | subj
## (Intercept) Residual
## StdDev: 139.2 179.4
##
## Fixed effects: rt ~ trial.type * mic
## Value Std.Error DF t-value p-value
## (Intercept) 757.6 42.23 12469 17.937 0.0000
## trial.typei 79.9 7.29 12469 10.957 0.0000
## micmicro 54.0 47.33 52 1.140 0.2594
## trial.typei:micmicro 17.8 8.17 12469 2.175 0.0297
## Correlation:
## (Intr) trl.ty micmcr
## trial.typei -0.068
## micmicro -0.892 0.061
## trial.typei:micmicro 0.061 -0.893 -0.068
##
## Standardized Within-Group Residuals:
## Min Q1 Med Q3 Max
## -4.6674 -0.5248 -0.1270 0.3272 11.4507
##
## Number of Observations: 12525
## Number of Groups: 54summary(m.het.mic.rea)## Linear mixed-effects model fit by REML
## Data: stroop.rt %>% filter(mic != "unknown")
## Subset: session == "rea"
## AIC BIC logLik
## 165469 165521 -82728
##
## Random effects:
## Formula: ~1 | subj
## (Intercept) Residual
## StdDev: 139.3 167.6
##
## Variance function:
## Structure: Different standard deviations per stratum
## Formula: ~1 | mic
## Parameter estimates:
## micro fomri
## 1.000 1.312
## Fixed effects: rt ~ trial.type + mic + trial.type:mic
## Value Std.Error DF t-value p-value
## (Intercept) 757.6 42.37 12469 17.881 0.0000
## trial.typei 79.9 8.94 12469 8.940 0.0000
## micmicro 54.0 47.44 52 1.137 0.2606
## trial.typei:micmicro 17.8 9.58 12469 1.856 0.0635
## Correlation:
## (Intr) trl.ty micmcr
## trial.typei -0.083
## micmicro -0.893 0.074
## trial.typei:micmicro 0.078 -0.934 -0.080
##
## Standardized Within-Group Residuals:
## Min Q1 Med Q3 Max
## -4.9965 -0.5329 -0.1253 0.3398 12.0000
##
## Number of Observations: 12525
## Number of Groups: 54summary(m.het.subj.rea)## Linear mixed-effects model fit by REML
## Data: stroop.rt %>% filter(mic != "unknown")
## Subset: session == "rea"
## AIC BIC logLik
## 160559 160997 -80220
##
## Random effects:
## Formula: ~1 | subj
## (Intercept) Residual
## StdDev: 138.2 174.8
##
## Variance function:
## Structure: Different standard deviations per stratum
## Formula: ~1 | subj
## Parameter estimates:
## 115825 123117 130114 130518 132017 138837
## 1.0000 2.5308 1.3131 0.4613 0.7545 1.1069
## 141422 150423 158136 165032 171330 178243
## 0.6145 0.8358 0.6950 0.7583 0.8991 1.4644
## 178647 178950 197449 203418 204319 233326
## 1.0869 0.9110 0.9354 0.6464 1.0887 0.8592
## 300618 317332 346945 393550 448347 580650
## 0.5252 1.0739 0.6571 0.7058 2.0377 1.7423
## 594156 601127 765864 814649 843151 849971
## 0.4918 1.1210 0.6392 1.0622 1.5554 1.1479
## 873968 877168 DMCC1328342 DMCC1596165 DMCC1971064 DMCC2442951
## 0.8275 1.5881 0.5940 1.7945 0.6483 0.3456
## DMCC2803654 DMCC2834766 DMCC3062542 DMCC4191255 DMCC4854075 DMCC5009144
## 0.6605 0.5432 0.8657 0.3361 0.3688 1.8097
## DMCC5195268 DMCC5775387 DMCC6371570 DMCC6418065 DMCC6671683 DMCC6705371
## 0.5287 1.0135 1.4468 0.6214 0.5124 0.4785
## DMCC6721369 DMCC8033964 DMCC8050964 DMCC9441378 DMCC9478705 DMCC9953810
## 0.4759 0.4892 0.7606 0.9560 0.8930 1.2308
## Fixed effects: rt ~ trial.type + mic + trial.type:mic
## Value Std.Error DF t-value p-value
## (Intercept) 760.5 41.95 12469 18.128 0.0000
## trial.typei 71.5 5.35 12469 13.364 0.0000
## micmicro 55.6 46.99 52 1.184 0.2417
## trial.typei:micmicro 12.5 5.86 12469 2.126 0.0336
## Correlation:
## (Intr) trl.ty micmcr
## trial.typei -0.050
## micmicro -0.893 0.045
## trial.typei:micmicro 0.046 -0.913 -0.049
##
## Standardized Within-Group Residuals:
## Min Q1 Med Q3 Max
## -6.2200 -0.6222 -0.1659 0.4498 10.6857
##
## Number of Observations: 12525
## Number of Groups: 54anova(m.hom.rea, m.het.mic.rea, m.het.subj.rea)## Model df AIC BIC logLik Test L.Ratio p-value
## m.hom.rea 1 6 165795 165839 -82891
## m.het.mic.rea 2 7 165469 165521 -82728 1 vs 2 327 <.0001
## m.het.subj.rea 3 59 160559 160997 -80220 2 vs 3 5014 <.0001plot(m.hom.pro)plot(m.het.subj.pro)d <- filter(stroop.rt, mic != "unknown", session == "pro")
plot(as.factor(d$trial.num), resid(m.het.subj.pro, type = "pearson"))plot(as.factor(d$run), resid(m.het.subj.pro, type = "pearson"))plot(as.factor(d$trial.type), resid(m.het.subj.pro, type = "pearson"))plot(as.factor(d$color), resid(m.het.subj.pro, type = "pearson"))plot(as.factor(d$word), resid(m.het.subj.pro, type = "pearson"))plot(as.factor(d$item), resid(m.het.subj.pro, type = "pearson"))plot(interaction(d$pc, d$trial.type), resid(m.het.subj.pro, type = "pearson"))plot(as.factor(d$iti), resid(m.het.subj.pro, type = "pearson"))plot(as.factor(d$mic), resid(m.het.subj.pro, type = "pearson"))plot(interaction(d$mic, d$trial.type, d$pc), resid(m.het.subj.pro, type = "pearson"))plot(
fitted(m.het.subj.pro), resid(m.het.subj.pro),
col = ifelse(d$mic == "micro", "black", ifelse(d$mic == "fomri", "firebrick", "grey50")),
pch = 16
)plot(
fitted(m.het.subj.pro), resid(m.het.subj.pro, type = "pearson"),
col = ifelse(d$mic == "micro", "black", ifelse(d$mic == "fomri", "firebrick", "grey50")),
pch = 16
)m.het.subj.trial.type.pro <- update(m.hom.pro, . ~ ., weights = varIdent(form = ∼ trial.type | subj))
summary(m.het.subj.trial.type.pro)## Linear mixed-effects model fit by REML
## Data: stroop.rt %>% filter(mic != "unknown")
## Subset: session == "pro"
## AIC BIC logLik
## 153344 153802 -76610
##
## Random effects:
## Formula: ~1 | subj
## (Intercept) Residual
## StdDev: 121.8 159.9
##
## Variance function:
## Structure: Different standard deviations per stratum
## Formula: ~trial.type | subj
## Parameter estimates:
## 115825 123117 130114 130518 132017 138837
## 1.0000 1.3684 0.9483 0.7232 1.1539 1.0273
## 141422 150423 158136 165032 171330 178243
## 1.0303 1.3350 0.6193 1.6752 0.9424 1.5864
## 178647 178950 197449 203418 204319 233326
## 1.1050 0.7494 1.1987 0.7306 1.8857 2.0556
## 300618 317332 346945 352738 393550 448347
## 0.5359 0.6395 0.8655 1.6850 0.7852 1.1969
## 580650 594156 601127 765864 814649 843151
## 2.0928 0.6314 1.1921 0.5608 0.9475 1.1984
## 849971 873968 877168 DMCC1328342 DMCC1596165 DMCC1971064
## 1.8200 0.7821 1.3117 0.4071 1.6184 0.8052
## DMCC2442951 DMCC2803654 DMCC2834766 DMCC3062542 DMCC4191255 DMCC4854075
## 0.4382 0.6706 0.7042 0.9748 0.6480 0.6854
## DMCC5009144 DMCC5195268 DMCC5775387 DMCC6371570 DMCC6418065 DMCC6661074
## 0.5089 0.5855 1.3942 0.9953 0.5380 0.7065
## DMCC6671683 DMCC6705371 DMCC6721369 DMCC7297690 DMCC8033964 DMCC8050964
## 0.5219 0.7398 0.5582 0.7739 0.5176 0.8544
## DMCC9441378 DMCC9478705 DMCC9953810
## 1.4022 0.8540 0.9206
## Fixed effects: rt ~ trial.type + mic + trial.type:mic
## Value Std.Error DF t-value p-value
## (Intercept) 735.5 35.44 11915 20.756 0.0000
## trial.typei 65.2 4.67 11915 13.940 0.0000
## micmicro 88.6 39.90 55 2.221 0.0305
## trial.typei:micmicro 8.7 5.45 11915 1.589 0.1121
## Correlation:
## (Intr) trl.ty micmcr
## trial.typei -0.087
## micmicro -0.888 0.078
## trial.typei:micmicro 0.075 -0.858 -0.091
##
## Standardized Within-Group Residuals:
## Min Q1 Med Q3 Max
## -4.9122 -0.6275 -0.1734 0.4266 10.1056
##
## Number of Observations: 11974
## Number of Groups: 57anova(m.het.subj.trial.type.pro, m.het.subj.pro)## Model df AIC BIC logLik
## m.het.subj.trial.type.pro 1 62 153344 153802 -76610
## m.het.subj.pro 2 62 153344 153802 -76610